Sequenced-based GWAS for linear classification traits in Belgian Blue beef cattle reveals new coding variants in genes regulating body size in mammals

Background Cohorts of individuals that have been genotyped and phenotyped for genomic selection programs offer the opportunity to better understand genetic variation associated with complex traits. Here, we performed an association study for traits related to body size and muscular development in intensively selected beef cattle. We leveraged multiple trait information to refine and interpret the significant associations. Results After a multiple-step genotype imputation to the sequence-level for 14,762 Belgian Blue beef (BBB) cows, we performed a genome-wide association study (GWAS) for 11 traits related to muscular development and body size. The 37 identified genome-wide significant quantitative trait loci (QTL) could be condensed in 11 unique QTL regions based on their position. Evidence for pleiotropic effects was found in most of these regions (e.g., correlated association signals, overlap between credible sets (CS) of candidate variants). Thus, we applied a multiple-trait approach to combine information from different traits to refine the CS. In several QTL regions, we identified strong candidate genes known to be related to growth and height in other species such as LCORL-NCAPG or CCND2. For some of these genes, relevant candidate variants were identified in the CS, including three new missense variants in EZH2, PAPPA2 and ADAM12, possibly two additional coding variants in LCORL, and candidate regulatory variants linked to CCND2 and ARMC12. Strikingly, four other QTL regions associated with dimension or muscular development traits were related to five (recessive) deleterious coding variants previously identified. Conclusions Our study further supports that a set of common genes controls body size across mammalian species. In particular, we added new genes to the list of those associated with height in both humans and cattle. We also identified new strong candidate causal variants in some of these genes, strengthening the evidence of their causality. Several breed-specific recessive deleterious variants were identified in our QTL regions, probably as a result of the extreme selection for muscular development in BBB cattle. Supplementary Information The online version contains supplementary material available at 10.1186/s12711-023-00857-4.

The results correspond to the MT-GWAS with traits related to buttock muscling (side view).The colors represent the LD level with the lead variant and the symbols indicate the predicted impact of the variant (• modifier, ♦ low impact, ▲moderate impact, ■ high impact).The positions of the genes are in the lower track.

Figure S2 .
Figure S2.Scatterplots for association levels for different traits for the QTL region on BTA4.The selected traits are those harboring a significant signal in the QTLR.Upper diagonal: scatterplots with p-values on a negative log10 scale.The color represents the LD level with the lead SNP (from the trait with the strongest association).Lower diagonal: scatterplots with signed t-values.The magenta circle denotes the lead variant.

Figure S3 .
Figure S3.Scatterplots for association levels for different traits for the QTL region on BTA5.The selected traits are those harboring a significant signal in the QTLR.Upper diagonal: scatterplots with p-values on a negative log10 scale.The color represents the LD level with the lead SNP (from the trait with the strongest association).Lower diagonal: scatterplots with signed t-values.The magenta circle denotes the lead variant.

Figure S4 .
Figure S4.Scatterplots for association levels for different traits for the QTL region on BTA6.The selected traits are those harboring a significant signal in the QTLR.Upper diagonal: scatterplots with p-values on a negative log10 scale.The color represents the LD level with the lead SNP (from the trait with the strongest association).Lower diagonal: scatterplots with signed t-values.The magenta circle denotes the lead variant.

Figure S5 .
Figure S5.Scatterplots for association levels for different traits for the QTL region on BTA16.The selected traits are those harboring a significant signal in the QTLR.Upper diagonal: scatterplots with p-values on a negative log10 scale.The color represents the LD level with the lead SNP (from the trait with the strongest association).Lower diagonal: scatterplots with signed t-values.The magenta circle denotes the lead variant.

Figure S6 .
Figure S6.Scatterplots for association levels for different traits for the QTL region on BTA18.The selected traits are those harboring a significant signal in the QTLR.Upper diagonal: scatterplots with p-values on a negative log10 scale.The color represents the LD level with the lead SNP (from the trait with the strongest association).Lower diagonal: scatterplots with signed t-values.The magenta circle denotes the lead variant.

Figure S7 .
Figure S7.Scatterplots for association levels for different traits for the QTL region on BTA19.The selected traits are those harboring a significant signal in the QTLR.Upper diagonal: scatterplots with p-values on a negative log10 scale.The color represents the LD level with the lead SNP (from the trait with the strongest association).Lower diagonal: scatterplots with signed t-values.The magenta circle denotes the lead variant.

Figure S8 .
Figure S8.Scatterplots for association levels for different traits for the QTL region on BTA23.The selected traits are those harboring a significant signal in the QTLR.Upper diagonal: scatterplots with p-values on a negative log10 scale.The color represents the LD level with the lead SNP (from the trait with the strongest association).Lower diagonal: scatterplots with signed t-values.The magenta circle denotes the lead variant.

Figure S9 .
Figure S9.Scatterplots for association levels for different traits for the QTL region on BTA26.The selected traits are those harboring a significant signal in the QTLR.Upper diagonal: scatterplots with p-values on a negative log10 scale.The color represents the LD level with the lead SNP (from the trait with the strongest association).Lower diagonal: scatterplots with signed t-values.The magenta circle denotes the lead variant.

Figure S10 .
Figure S10.Regional association plot for the QTLR associated to recessive genetic defects.The colors represent the LD level with the lead variant and the symbols indicate the predicted impact of the variant (• modifier, ♦ low impact, ▲moderate impact, ■ high impact).The positions of the genes are in the lower track.Upper panel: mapping on BTA3 (MT-GWAS on muscular development traits); middle panel: mapping on BTA19 (MT-GWAS on body size traits); lower panel: mapping on BTA25 (MT-GWAS on muscular development traits).

Figure S11 .
Figure S11.Regional association plot for the QTLR on BTA26.The colors represent the LD level with the lead variant and the symbols indicate the predicted impact of the variant (• modifier, ♦ low impact, ▲moderate impact, ■ high impact).The positions of the genes are in the lower track.Upper panel: ST-GWAS for height; middle panel: MT-GWAS on muscular development traits; lower panel: MT-GWAS on body size.

Figure S12 .
Figure S12.Regional association plot for the QTLR on BTA6.The results correspond to the MT-GWAS with traits related to body size.The colors represent the LD level with the lead variant and the symbols indicate the predicted impact of the variant (• modifier, ♦ low impact, ▲moderate impact, ■ high impact).The positions of the genes are in the lower track.

Figure S13 .
Figure S13.Regional association plot for the QTLR on BTA5.The results correspond to the MT-GWAS with traits related to body size.The colors represent the LD level with the lead variant and the symbols indicate the predicted impact of the variant (• modifier, ♦ low impact, ▲moderate impact, ■ high impact).The positions of the genes are in the lower track.

Figure S14 .
Figure S14.Regional association plot for the QTLR on BTA18.The results correspond to the MT-GWAS with traits related to body size.The colors represent the LD level with the lead variant and the symbols indicate the predicted impact of the variant (• modifier, ♦ low impact, ▲moderate impact, ■ high impact).The positions of the genes are in the lower track.

Figure S15 .
Figure S15.Regional association plot for the QTLR on BTA23.The results correspond to the MT-GWAS with traits related to body size.The colors represent the LD level with the lead variant and the symbols indicate the predicted impact of the variant (• modifier, ♦ low impact, ▲moderate impact, ■ high impact).The positions of the genes are in the lower track.

Figure S16 .
Figure S16.Regional association plot for the conditional mapping in QTLR on BTA3, BTA14, BTA16 and BTA26.The left panels represent the initial GWAS whereas the right panels correspond to conditional GWAS in which the candidate variants are fitted as covariate.The colors represent the LD level with the lead variant.The positions of the genes are in the lower track.a) GWAS for top muscling on BTA3, b) GWAS for pelvis width on BTA14, c) GWAS for shoulder muscling on BTA16, and d) GWAS for pelvis length on BTA26.

Figure S17 .
Figure S17.Regional association plot for the conditional mapping in QTLR on BTA5 and BTA23.The left panels represent the initial GWAS whereas the right panels correspond to conditional GWAS in which the candidate variants are fitted as covariate.The colors represent the LD level with the lead variant.The positions of the genes are in the lower track.a) GWAS for pelvis length on BTA5, b) GWAS for length on BTA5, c) GWAS for height on BTA23, and d) GWAS for pelvis length on BTA23.

Figure S18 .
Figure S18.Regional association plot for the conditional mapping in QTLR on BTA19.The left panels represent the initial GWAS whereas the right panels correspond to conditional GWAS in which the candidate variants are fitted as covariate.The colors represent the LD level with the lead variant.The positions of the genes are in the lower track.GWAS for a) length, b) top muscling, c) rump, and d) height.

Figure S19 .
Figure S19.Regional association plot for the conditional mapping in QTLR on BTA25.The left panels represent the initial GWAS whereas the right panels correspond to conditional GWAS in which the candidate variants are fitted as covariate.The colors represent the LD level with the lead variant.The positions of the genes are in the lower track.GWAS for a) rump, and b) buttock muscling (side view).The boxes indicate the two regions achieving similar association levels and included in the credible set.

Figure S20 .
Figure S20.Regional association plot for the conditional mapping for the QTLR on BTA25.The results correspond to the MT-GWAS with traits related to buttock muscling (side view).The colors represent the LD level with the lead variant and the symbols indicate the predicted impact of the variant (• modifier, ♦ low impact, ▲moderate impact, ■ high impact).The positions of the genes are in the lower track.